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Abstract. We show that the MHD theory that explains the solar tachocline by an effect of the magnetic field can work with the 
decay modes of a fossil field in the solar interior if the meridional flow of the convection zone penetrates slightly the radiative 
zone beneath. An equatorward flow of about 10 m/s penetrating to a maximum depth of 1000 km below the convection zone 
is able to generate almost horizontal field lines in the tachocline region so that the internal field is almost totally confined to 
the radiative zone. The theory of differential solar rotation indeed provides meridional flows of about 10 m/s and a penetration 
depth of < 1000 km for viscosity values that are characteristic of a stable tachocline. 
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1. Introduction 

While rotation of the solar radiative core is almost uniform, 
the convection zone rotates differentially (Wilson et al. 1997 
Kosovichev et al. 1997 Shou et al.[l998 1. The transition from 
inhomogeneous to rigid rotation occurs in a thin almost spheri- 
cal layer named 'tachocline', whose parameters are well known 
from helioseismology. The tachocline thickness is estimated 
to be < O.O5/? , its midpoint radius is (0.692 ± 0.005) R Q , and 
it is slightly prolate in shape (Kosovichev 119961 Antia et al. 
1998 Charbonneau et al. 1999b I. The tachocline lies totally, 
or at least mostly, beneath the base of convection zone at 
R m = O.713fl (Christensen-Dalsgaard et al. AWTl Basu & 



Antia 1997), at the very top of the radiative core. 

The existence of tachocline is evidence of a link between 
low and high latitudes that should be present just beneath the 
convection zone in order to suppress the latitudinal differential 
rotation. Similarly, uniform rotation in the radius beneath the 
tachocline indicates that a coupling in radius should be present 
there. Otherwise, rotational braking of the Sun would leave a 
rapidly rotating core. Among several theoretical concepts sug- 
gested for the solar tachocline, only one that explains it by the 
effect of a weak internal magnetic field can provide both links 
simultaneously (Charbonneau & MacGregor 1993 MacGregor 
& Charbonneau 1999). The uniform rotation of the core was 
even considered as evidence of an internal field. Another argu- 
ment by Cowling ( 1945) is that the time of resistive decay for 
the hypothetical internal field exceeds the solar age. 



It has been confirmed by numerical modeling that even 
a weak (~ 10~ 3 G) poloidal field can produce the tachocline 
(Rudiger & Kitchatinov 119971 MacGregor & Charbonneau 
1999). However, this is the case only if the field geometry sat- 
isfies certain constraint. Ferraro's dl937t law on the constancy 
of angular velocity along the field lines requires the field to be 
almost horizontal inside the tachocline region (except, perhaps, 
for the mid latitudes where angular velocity at the base of con- 
vection zone equals that of the bulk of the core). 

Why the field should be confined inside the radiative core 
has remained, however, uncertain. The main idea of the present 
paper is that the required field geometry may result from the 
influence of the meridional flow penetrating the radiative core 
from the convection zone. If the electric conductivity in the ra- 
diative core is much higher than in the convection zone, then 
(and only then) the flow has massive electrodynamical conse- 
quences, and the poloidal magnetic field becomes parallel to 
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the meridional flow (cf. Mestel [1999l ). 

The flow - produced in the convection zone by centrifugal 
and barocline forcing - only slightly penetrates the radiative 
core (Gilman & Miesch 2004; Kitchatinov & Rudiger 12005V 
Its induction can thus be included via the boundary condi- 
tion. In this paper we formulate the boundary condition for a 
poloidal field on the top of the radiative core with account for 
the shallow penetration of meridional flow. With such a modi- 
fied boundary condition, the modes of the internal poloidal field 
with the longest decay times are computed. The long-living 
modes show the required confined structure. The thus defined 



2 



L. L. Kitchatinov and G. Riidiger: Magnetic solar tachocline 



internal field is then used in tachocline computations to confirm 
that a slender tachocline can indeed be found. 



2. Penetration 

The penetration of the meridional flow beneath the solar con- 
vection zone was computed by Kitchatinov & Riidiger (2005 1. 
The main results are shown in Fig.^ The meridional velocity 
in the penetration region shows reversals with increasing depth, 
where the flow amplitude reduces largely after each reversal. 
The penetration depth D pen in Fig. [Qhas been defined as the 
distance from the base of the convection zone to the location of 
the first reversal. 
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Fig. 1. Left: penetration depth as a function of the core viscos- 
ity. Right: meridional velocity profile in the penetration region 
at 45° latitude for v in = 1.1 ■ 10 9 cm s . Negative values mean 
equatorward flow. The vertical dotted line marks the bottom of 
the convection zone. 



The penetration results from viscous drag imposed by the 
meridional flow at the base of the convection zone on the fluid 
beneath, and this is opposed by the Coriolis force. The Ekman 



balance results in the estimation D 



pen 



Vvin/Q (Gilman & 
Miesch 2 0041 v m is the viscosity in the penetration region, Q. its 
basic angular velocity). The plot of Fig.^Oeft) is approximated 
well by Dpen - 2.3 -\Jv m /Q or, equivalently, by 



£> P en = 1.4Vvta' 10 3 cm 



(1) 



for solar parameters. Thus, the penetration distance is about 
100 m for microscopic viscosity and remains smaller than 
1 Mm for any reasonable value of eddy viscosity. The penetra- 
tion depth is small even compared to the tachocline thickness, 

^pen ^ Aach- 

This shallow penetration may still strongly influence the 
geometry of an internal poloidal field. The penetration makes 
a shear flow on the top of the radiative core. The time r s 
D p /u m of latitudinal field shearing from the radial one is small 
compared to the time t&b = D 2 en /j]i n of diffusive escape of the 
field from the penetration region. The time ratio, t&b/Ts, gives 
the magnetic Reynolds number, 



(2) 



where Pm is the magnetic Prandtl number. The value u m =s 
10 ms -1 is used for the meridional velocity. With microscopic 



diffusivities, the number is large, Rm ~ 10 3 , and it remains 
higher than this value with eddy diffusivities up to 10 6 cm 2 s _1 . 
It sinks to unity for eddy diffusivities of about 10 12 cm 2 s _1 , so 
that only in this case does the influence of a penetrating merid- 
ional flow on the magnetic field geometry turn weak. 

The shearing of the poloidal field is probably the only pro- 
cess where this penetration is important. The advection time, 
T a dv = R'm/u m , is much longer than the diffusion time Tdiff, 

D n 



Tdiff 
^"adv 



'pen 



Rm 0.03 • Pm <sc 1, 



(3) 



independent of whether microscopic or eddy diffusion applies. 
Riidiger et al. ( 2005 1 argued that the tachocline region should 
be stable in the hydrodynamic regime so that microscopic dif- 
fusivities have to be used. In any case the diffusion time is so 
short that the penetration layer is not dynamo-relevant (Riidiger 
et al. 2005). 

3. The model 

Axial symmetry for both magnetic and velocity fields is as- 
sumed. In spherical coordinates (r, 9, <p), the fields can be writ- 
ten as 

1 dtfi -1 dip 



B 



pr 2 sin 9 39' pr sin 9 dr 
1 dA -1 dA 



r sin 9 Q 



r 2 sin 9 89 r sin 9 dr 



B 



(4) 
(5) 



in terms of the stream function (t//) of the meridional flow and 
the potential (A) of the poloidal field. Here Q is the angular 
velocity and B the toroidal field. 
The poloidal field equation, 

dA _ ugdA dA d 2 A T] . „<9 / 1 dA] 
~dt ~ 



d9 



dr dr L 



77 d [ 

sm9— — 

r 2 d9 \ sin 9 d9 , 



(6) 



decouples from other equations of the problem. As the eddy 
diffusivity in the convection zone is large compared to rj m , the 
vacuum boundary condition can be applied to the top of the 
core. The (nonlocal) condition, 



dA _ I8A\ 
3r~\drL 



at r = R in , 



(7) 



is usually expressed in terms of the Legendre polynomial ex- 
pansion 



A(r,9) = sinO } j A n {r)P l n (cos 0), 
dA\ 



n=\ 

sin# 



oa \ sin w . , n i , 

= \nA„{r)P\{cos9). 

/vac h— 1 



(8) 



The internal magnetic field can be defined by solving 
Eq. with the boundary condition Q. The meridional flow is, 
however, only present in the very thin penetration layer (Fig.Q. 
Instead of resolving this layer explicitly, we integrate Eq. 
across the layer to obtain a new boundary condition for the 
poloidal field for the base of the penetration layer that is still 
above the tachocline. The effect of penetration is then included 
in the reformulated boundary condition. 
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3.1. Poloidal field boundary condition 

Inside the penetration layer, the advection time (~ 10 yr) is 
short compared to the latitudinal diffusion time, R 2 /j] m ~ 
10 9 yr. The left part of Eq. (|6j and its last term can thus be 
neglected compared to the first term on the right. 

Also, the radial velocity in the penetration layer is small 
compared to the latitudinal velocity, u r ~ ueD pen /Ri n (Gilman 
& Miesch 2004 1. By neglecting these small terms the equation 
for the penetration region reads 



1 



diffdA 



8 2 A 



dr 2 rjpR 2 sin 6» dr 80 



0, 



(9) 



where ue is expressed in terms of the stream function of Eq. (0}. 
The stream function varies from zero at the base of the layer 
to a finite value at its top. Divergence-free of magnetic field 
requires a variation in the radial field, B r - (dA/dff)/{R 2 n sin 9), 
across the layer be small. Integration of Eq. © across the layer 
then yields 



8A 
dr 



HRin,0) OA 
T] in pR 2 sin 80 



8A\ 



(10) 



where if/{R m , 0) is the stream function profile at the base of the 
convection zone. The stream function scales as 



(ID 



where \j/ is a dimensionless function of order one. Substitution 
of dl It into (II Oi leads to the boundary condition sought for the 
poloidal field, 



D 8A £(0)dA 
Rm~Z Rm — 

dr sm0 00 



(dA\ 



(12) 



which accounts for the shallow penetration of the meridional 
flow from the convection zone into the radiative interior. 



3.2. Model equations 

Any meridional flow in the bulk of radiative core is neglected. 
The characteristic time of the Eddington-Sweet circulation ex- 
ceeds the solar age (Tassoul 2000i. The slow circulation may 
affect the internal field structure but the effect remains small 
(Garaud l2"UU2l . 

The global modes of the internal poloidal field of the core 
can be computed by solving the eigenvalue equation, 



A 

T 



(Pa 

' dr 2 



I^TJ + — sin#— I 



i dA 



80 \sin0 00 



(13) 



with the boundary condition d!2i . The other condition requires 
the solution to be regular at the center. We are interested in the 
solutions with the longest decay times, t. 

The radial profiles of the microscopic diffusivities and den- 
sity were defined after the solar structure model by Stix & 
Skaley (119901 see Riidiger & Kitchatinov 1 19971 . All the com- 
putations concern the simplest example of a stream function, 



ij/(6) = -s\n0P\ (cos 0), 



(14) 



where P x 2 is the normalized Legendre polynomial. 
Equation ( I14i describes a flow at the base of the convec- 
tion zone from the poles to the equator, which appeared in all 
our previous models for the solar differential rotation. 

To describe the poloidal field geometry, the 'escape param- 
eter' , 



max|A (r, 0) \ r=Riri 
max|A (r, 0) \ r < R ' 



(15) 



is used. The value of 27rA (r, 0) equals the magnetic flux per- 
vading a surface bounded by the longitudinal circle of constant 
r and 0. The <50-parameter 11151 estimates the ratio of the mag- 
netic flux through the surface of the core to the characteristic 
value of the flux inside the core. The poloidal fields with small 
5(f> are expected to produce the tachocline. Linear equation dl 31 
leaves the amplitude of the field indefinite. The amplitude, Bo, 
defined as the maximum strength of poloidal field within the 
core, remains a free parameter of the model. 

With a given poloidal field, the tachocline can be computed 
by solving the equation system for the toroidal field and angular 
velocity, 

7] d I 1 d(Bsin0)\ d I d(Br)\_ 8D.8A 8Q.8A 
7~80\l\n0 80 j + Yr y dr j ~ ~80~~8r~ ~ ~dr"d~0'' 



pv 8 , 5Q\ 13, ffl 
At-— sin 3 0— + -r— \r 4 pv—\ 
sin 3 80 \ 80 r 2 8r* > ' 1 



1 



47ir 2 sin 3 \ dr 



8A8(B sinfl) 

r— sinfl 

80 



dr ) 
8A8(Br) 



80 dr 



(16) 



The neglect of the meridional flow in these equations is justi- 
fied by the small value of parameter Q. The upper boundary 
conditions require a zero toroidal field and the angular velocity 
profile 



Q. = 2.9(l -0.15 cos 2 9) ^rads" 



(17) 



favored by the helioseismology for the bottom of the convec- 
tion zone (Charbonneau et al. 1999a). 

We applied a Legendre polynomial expansion in to Eqs. 
dl 3i and d!6l > and solved the resulting system of ordinary dif- 
ferential equations numerically with the relaxation method. 

4. Results and discussion 

Figures 13 and |3 illustrate the effect of the penetrating merid- 
ional flow on the structure of the longest-living dipolar mode 
of an internal poloidal field. The internal field of Fig.|5]com- 
puted without penetration (Rm = 0) has an 'open' structure. 
The field computed with Rm = 1000 has the 'confined' geom- 
etry required for tachocline formation. Figure [3] shows the de- 
pendence of the escape parameter d!5l on Rm. For Rm = 1000, 
less than 1 % of magnetic flux of the dipolar eigenmode escapes 
the core. 

The decay time of the longest-living mode of internal field 
exceeds the solar age. The lifetime decreases with Rm but ap- 
proaches a constant value of about 5.4 Gyr for large Rm. The 
characteristic time of the tachocline formation, ~ R\ n /VA ~ 
10 4 /Bq yr (Va is Alfven velocity, Bo is the field amplitude 
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Fig. 2. Field lines of the longest-living modes of poloidal field 
without (Rm = 0, left) and with penetration (Rm = 1000, 
right). The outer circle is the base of the convection zone. 
Lifetimes of the fields in Gyr are given at the top. 
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Fig. 3. Geometry parameter (1151 for the longest-living dipo- 
lar mode of the internal field as a function of the magnetic 
Reynolds number (0. The field becomes increasingly confined 
to the core as Rm grows. 



in Gauss) is much shorter. This justifies using steady poloidal 
fields in the tachocline computations. 

Figure|4]shows the angular velocity distributions inside the 
core computed for poloidal fields of Fig. [2 As expected, only 
the confined field produces a tachocline. The dipolar eigen- 
mode for the solar value of Rm = 1000 is 'sufficiently con- 
fined' to do so. Already Rm = 100 suffices for the tachocline 
formation. We conclude that the shallow penetration of the 
meridional flow into the radiative core influences the inter- 
nal field geometry so strongly that it becomes appropriate for 
tachocline formation. We suggest that a layer of less than 1000 
km beneath the convection zone, where a meridional flow on 
the order of 1-10 m/s enters the radiative core, is responsible 
for the existence of the impressive phenomenon of the solar 
tachocline (which implies that the bulk of the solar core rotates 
rigidly). 

Analytical estimations suggest that the (fractional) 
tachocline thickness is controlled solely by the Hartmann 
number, Ha = V A Rm/ >/W^> 




Fig. 4. Angular velocity isocontours inside the radiative zone 
computed with the internal fields of Fig. 13 Left: Rm = 0, right: 
Rm = 1000. B = 3 mG. 



(Riidiger & Kitchatinov FT9971 . After Eq. (|T8}, even a weak 
internal field can produce the tachocline. The angular veloc- 
ity distributions of Fig. 0] were computed for a field ampli- 
tude of 3 mG. Figure |5] shows the fractional depth inside the 
core at which the angular velocity difference between equa- 
tor and poles drops e-times. When the field is not too weak, 

— 1/2 

the tachocline thickness varies as Aach ~ B {) in agreement 
with (11 8t . The thickness falls below 4% of the solar radius at 
Bo = 10~ 2 G. The toroidal field in the present model is almost 
independent of Bo, and its amplitude is not large, B t01 ~ 100 G. 
The maximum toroidal field within the core is about 140 G for 
Rm = 1000. 
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-1/2 



(18) 



Fig. 5. Dependence of the thickness of the tachocline on the 
amplitude, Bo, of the internal poloidal field. The thickness is 
defined as the depth of an exponential reduction in the pole- 
equator difference of the angular velocity. 



The dipolar mode of Fig. |2] has the longest lifetime. The 
following modes in order of decreasing lifetime are shown 
in Fig. [6] The lifetimes are shorter than the solar age but 
not by much. The modes may be mixed in the internal field. 
Tachocline computations were made for all modes of Fig. [6] 
and all the field configurations lead to a slender tachocline. 

It is important for this model that the magnetic Reynolds 
number be large, which is true for microscopic diffusiv- 
ities or for eddy diffusivities that are not too large, up to 
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Fig. 6. Eigenmodes of the internal field following the dipolar 
mode of Fig.[5]in order of decreasing lifetime. The lifetimes in 
Gyr are given at the top. Rm = 1000. The geometry parameter 
IT3I varies as 0.0051, 0.0073, 0.0026 from left to right. 



about 10 8 cm 2 s . The tachocline should thus be stable or only 
mildly turbulent in order to allow sufficiently large Rm. 

Another necessary condition is the presence of a meridional 
flow of about 10 ms"' at the bottom of the convection zone. 
This flow should, of course, be present in the deep convection 
zone in order to penetrate below its base. Such a flow is also a 
key ingredient in the advection-dominated solar dynamo mod- 
els (Wang et al. lT99Tl Choudhuri et al. lT995l Dikpati & Gilman 
2001 Bonanno et al. 2002 1. The flow was predicted by theoret- 
ical modeling (Kitchatinov & Riidiger 1999 Miesh et al. 2000 1, 
but not yet confirmed by any other means. Probing for the deep 
meridional flow may be a challenge for helioseismology. 

The polar cusp in the tachocline of Fig. 0] (right) is a con- 
sequence of the assumed axial symmetry. This symmetry is, of 
course, an idealization. We expect, though so far only on qual- 
itative grounds, that i) the meridional flow penetration can also 
confine a nonaxisymmetric field to the core and ii) the con- 
fined nonaxisymmetric field can also produce the characteristic 
tachocline structure. An important question is still open as to 
whether the axisymmetric solutions are stable against nonax- 
isymmetric disturbances. Purely toroidal fields (Tayler 1973 1, 
as well as the poloidal fields (Wright [TW3J Markey & Tayler 
I1973> . are known to be unstable near their neutral lines if the 
fields are sufficiently strong, so that the characteristic time of 
an Alfven wave passage across the star is shorter than the rota- 
tion period. The weak fields of the present model belong to the 
opposite limit. 
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